The aim of this project is to develop a generic dashboard, e.g. a R shiny application to 1) inspect model inputs and outputs, 2) visualise the original inputs and outputs, 3) investigate the relationship between model inputs and outputs through metamodelling and data visualisation methods, 4) save the performed analyses.
Package with supporting functions can be found here: https://github.com/Xa4P/pacheck.
The envisoned R Shiny app will have the following tabs. I’ve developed some functions already to implement. Currently, I’ve only considered 2 strategies (“intervention” (_in) and “comparator” (_comp)). What is shown in this document is based on the probabilistic analysis of a (toy) 3-states Health state transition model (Progression-free (PF), Progressed disease(PD), Dead(D)) depicted in Figure 1. The intervention is only having an effect on the probability of progression, and incurs costs in the Progression-Free health state.
Content:
df_pa <- calculate_nb(df = df_pa,
e_int = "QALY_int",
e_comp = "QALY_comp",
c_int = "Costs_int",
c_comp = "Costs_comp",
wtp = 80000)# calculate net benefits
To examine whether cost inputs are always positive for instance
df <- generate_sum_stats(df_pa)
kable(df)
| Parameter | Mean | SD | Percentile_2.5th | Percentile_97.5th | Minimum | Maximum |
|---|---|---|---|---|---|---|
| p_pfspd | 0.150 | 0.035 | 0.088 | 0.226 | 0.049 | 0.316 |
| p_pfsd | 0.100 | 0.030 | 0.049 | 0.165 | 0.023 | 0.244 |
| p_pdd | 0.201 | 0.040 | 0.129 | 0.283 | 0.073 | 0.366 |
| p_dd | 1.000 | 0.000 | 1.000 | 1.000 | 1.000 | 1.000 |
| rr | 0.753 | 0.068 | 0.630 | 0.895 | 0.535 | 1.033 |
| u_pfs | 0.816 | 0.065 | 0.675 | 0.924 | 0.491 | 0.974 |
| u_pd | 0.587 | 0.103 | 0.382 | 0.776 | 0.207 | 0.898 |
| u_d | 0.000 | 0.000 | 0.000 | 0.000 | 0.000 | 0.000 |
| c_pfs | 998.950 | 202.516 | 603.870 | 1392.049 | 202.083 | 1771.015 |
| c_pd | 1993.689 | 401.306 | 1203.104 | 2782.842 | 388.754 | 3707.492 |
| c_d | 0.000 | 0.000 | 0.000 | 0.000 | 0.000 | 0.000 |
| c_thx | 10000.807 | 98.670 | 9807.326 | 10194.863 | 9632.935 | 10341.933 |
| QALY_comp | 4.005 | 0.743 | 2.727 | 5.644 | 1.939 | 7.889 |
| QALY_int | 4.309 | 0.824 | 2.894 | 6.105 | 2.104 | 8.433 |
| Costs_comp | 7223.596 | 1593.274 | 4473.077 | 10703.240 | 2425.117 | 17054.204 |
| Costs_int | 38841.599 | 7233.474 | 26434.052 | 54575.070 | 19181.834 | 78244.121 |
| Inc_QALY | 0.304 | 0.170 | 0.050 | 0.703 | -0.291 | 1.562 |
| Inc_Costs | 31618.003 | 6576.688 | 20379.534 | 46057.098 | 15154.291 | 66447.151 |
| NMB_int | 305872.059 | 60085.800 | 202019.061 | 436642.750 | 148476.152 | 609534.328 |
| NMB_comp | 313195.478 | 58530.830 | 212540.838 | 442388.664 | 151753.123 | 614042.727 |
| iNMB | -7323.419 | 10571.658 | -24480.486 | 16838.699 | -52000.775 | 61320.703 |
| NHB_int | 3.823 | 0.751 | 2.525 | 5.458 | 1.856 | 7.619 |
| NHB_comp | 3.915 | 0.732 | 2.657 | 5.530 | 1.897 | 7.676 |
| iNHB | -0.092 | 0.132 | -0.306 | 0.210 | -0.650 | 0.767 |
rm(df)
To visually investigate the parameter distributions
p_1 <- vis_1_param(df = df_pa,
param = "u_pfs",
binwidth = NULL,
type = "histogram",
dist = c("norm", "beta", "gamma", "lnorm"))
## Loading required package: fitdistrplus
## Loading required package: MASS
##
## Attaching package: 'MASS'
## The following object is masked from 'package:plotly':
##
## select
## Loading required package: survival
p_1
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
p_2p <- vis_2_params(df = df_pa,
param_1 = "u_pfs",
param_2 = "u_pd",
slope = 1,
check = "param_2 > param_1")
## [1] "P(TRUE): 3 %"
p_2p
To cross check with parameters reported in documentation/report/ journal article, as an implementation check
p_2 <- vis_1_param(df = df_pa,
param = "u_pfs",
binwidth = NULL,
type = "density",
dist = c("norm", "beta", "gamma", "lnorm"))
p_2
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
Questions:
p_3 <- plot_ice(df = df_pa,
param_1 = "Inc_QALY",
param_2 = "Inc_Costs",
wtp = 80000)
## Loading required package: scales
p_3
p_3_interact <- ggplotly(p_3)
p_3_interact
m_ceac <- calculate_ceac(df = df_pa,
e_int = "QALY_int",
e_comp = "QALY_comp",
c_int = "Costs_int",
c_comp = "Costs_comp")
df_ceac <- as.data.frame(m_ceac)
plot_ceac(df = df_ceac,
wtp = "WTP_threshold")
## Loading required package: reshape2
Can use the function above!
###d. Convergence graph of outcomes
plot_convergence(df = df_pa,
outcome = "iNMB"
)
lm_rr <- fit_lm_metamodel(df = df_pa,
x = "rr",
y = "iNMB")
lm_pred <- unlist(predict(lm_rr, data.frame(rr = df_pa$rr)))
df_obs_pred <- data.frame(
Values = df_pa$rr,
Observed = df_pa$iNMB,
Predicted = lm_pred
)
ggplot(data = df_obs_pred, aes(x = Values, y = Observed)) +
geom_point(shape = 1, colour = "lightgrey") +
geom_smooth(method = "lm") +
theme_bw()
## `geom_smooth()` using formula 'y ~ x'
plot(lm_rr)
lm_full <- fit_lm_metamodel(df = df_pa,
x = c("rr", "u_pfs", "u_pd", "c_pfs", "c_pd", "c_thx", "p_pfspd", "p_pfsd", "p_pdd"),
y = "iNMB")
lm_pred_full <- predict(lm_full, data.frame(rr = df_pa$rr,
u_pfs = mean(df_pa$u_pfs),
u_pd = mean(df_pa$u_pd),
c_pfs = mean(df_pa$c_pfs),
c_pd = mean(df_pa$c_pd),
c_thx = mean(df_pa$c_thx),
p_pfspd = mean(df_pa$p_pfspd),
p_pfsd = mean(df_pa$p_pfsd),
p_pdd = mean(df_pa$p_pdd)))
df_obs_pred_full <- data.frame(
Values = df_pa$rr,
Observed = df_pa$iNMB,
Predicted = lm_pred_full
)
ggplot(data = df_obs_pred_full, aes(x = Values, y = Observed)) +
geom_point(shape = 1, colour = "lightgrey") +
geom_line(data = df_obs_pred_full, aes(x = Values, y = Predicted), colour = "blue") +
theme_bw()
plot(lm_full)
df_dsa <- dsa_lm_metamodel(df = df_pa,
lm_metamodel = lm_full)
plot_tornado(df = df_dsa,
df_basecase = df_pa,
outcome = "iNMB")
## Loading required package: tidyverse
## -- Attaching packages --------------------------------------- tidyverse 1.3.0 --
## v tibble 3.0.5 v dplyr 1.0.3
## v tidyr 1.1.2 v stringr 1.4.0
## v readr 1.4.0 v forcats 0.5.1
## v purrr 0.3.4
## -- Conflicts ------------------------------------------ tidyverse_conflicts() --
## x readr::col_factor() masks scales::col_factor()
## x purrr::discard() masks scales::discard()
## x dplyr::filter() masks plotly::filter(), stats::filter()
## x dplyr::lag() masks stats::lag()
## x dplyr::select() masks MASS::select(), plotly::select()
testthat? (automatic check of costs are positive for instance, or utility values between 0-1).